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Abstract 

A nonlinear internal model control (NIMC) strategy that 
incorporates the nonlinear model structure and the estimator 
dynamics in the control law is presented for the control of 
complex dynamical systems characterized by input-output 
multiplicities, nonlinear oscillations and chaos. A model 
based estimator is designed to provide the unmeasured 
process states that capture the fast changing nonlinear 
dynamics of the process to incorporate in the controller. The 
estimator uses the mathematical model of the process in 
conjunction with the known process measurements to 
estimate the states. The design and implementation of the 
estimator supported NIMC strategy is studied by choosing 
two typical continuous non-isothermal nonlinear processes, 
a chemical reactor and a polymerization reactor, which show 
rich dynamical behaviour ranging from stable situations to 
chaos. The results evaluated under different conditions show 
the superior performance of the estimator based NIMC 
strategy over the conventional controllers for the control of 
complex nonlinear processes. 
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Introduction 

Control of nonlinear systems exhibiting complex 
dynamic behaviour is a challenging task because such 
systems present a variety of behavioural patterns 
depending on the values of their physical parameters 
and intrinsic features. Among these systems, 
continuous chemical reactors and continuous 
polymerization reactors have received special 
attention due to their complicated dynamics and 
economic importance. Depending on the parameter 
values, these systems can be operated at steady state 
or present oscillatory and chaotic motions. Several 
researchers have analysed continuous stirred tank 
chemical reactors to characterize the phenomena of 
multiple steady states, simple and modulated 



oscillations and chaos [Uppal, Ray and Poore, 1976, 
Schmitzet al., 1979, Doedel and Heinemann, 1983, 
Blanco and Bandoni, 2007]. Certain continuous 
polymerization reactors exhibit highly nonlinear 
dynamics because of the complicated reaction 
mechanisms associated with the large number of 
interactive reactions. Many researchers have 
investigated the existence of steady state multiplicities; 
self-sustained oscillations and chaos in continuous 
solution polymerization reactors [Jaisinghani and Ray, 
1977, Schmidtet al., 1984, Teymour and Ray, 1989, 
Flores-Tlacuahuac et al., 2005]. 

The nonlinear dynamical systems due to variations in 
the operating parameter space can lead to exhibit 
open-loop unstable behaviour and some chemical 
systems are more sensitive to these variations where 
small perturbation in the operating conditions can 
cause the stable operation to degenerate into 
oscillatory and chaotic motion. The unstable, 
oscillatory and chaotic phenomenon displayed by the 
chemically reacting systems has desirable as well as 
undesirable features. The desirable feature of multiple 
steady state condition is that one of the unstable 
steady state may correspond to higher rate of 
reaction/yield and this becomes the preferred state to 
be operated on it for enhancing the process 
performance. The desirable feature of chaos is that it 
enhances mixing and chemical reactions and provides 
a vibrant mechanism for transport of heat and mass. 
On the other hand, the intrinsic features of the reacting 
systems with the interactive influence of chemical or 
thermal energy may cause irregular dynamic 
behaviour leading to degraded performance. In such 
situations, chaotic behaviour is considered as 
undesirable and should be avoided. Chaotic processes 
show extreme sensitivity to initial conditions and the 
process trajectories can diverge exponentially making 
the task of control more difficult. However, chaos 
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offers great flexibility to operate chemical systems 
because there are an infinite number of unstable 
periodic orbits (UPO) embedded in a chaotic attractor, 
in which one can choose a specific UPO along with its 
time averaged pre-specified performance defined in 
terms of conversion, yield or selectivity and the 
process can be stabilized to operate it on the chosen 
UPO. Various conventional and advanced techniques 
have been reported for the control of open-loop 
unstable chemical and polymerization reactors and the 
reactors that show oscillatory and chaotic dynamics. 
As far as unstable reactors are concerned, Choi and 
Ray [1988] and Chan et al. [1993] have employed 
standard PID controllers to stabilize open-loop 
unstable steady state conditions of an olefin 
polymerization reactor, however, Takamatsu et al. 
[1986] and Kwag and Choi [1994] showed that such 
controllers may also pose serious control problems if 
the polymerization reactor presents significant 
nonlinear behaviour. Dittamar et al. [1991] reported 
that linear predictive controllers might lead to proper 
control of open-loop unstable polymerization reactor. 
Various nonlinear model predictive control strategies 
have been successfully applied to the control of open- 
loop unstable CSTR systems [Hidalgo and Brosilow, 
1990, Hapoglu, 2002, Venkateswarlu and Venkat Rao, 
2005]. Nonlinear processes characterized with chaotic 
dynamic behaviour can be controlled either by 
stabilizing the UPOs in the neighbourhood of the 
desired UPO or stabilizing the dynamic trajectory of 
the system exactly at the unstable steady state. Various 
controllers based on ad-hoc, conventional and 
advanced approaches have been reported to stabilize 
the chaotic systems. In the ad-hoc approach, the 
properties of chaos are used to develop ad-hoc 
methods without considering conventional control 
techniques. The most commonly used method is the 
OGY method by which one of the unstable orbits in 
the chaotic attractor can be stabilized via small time 
dependent perturbations to a system parameter [Ottet 
al., 1990, So and Ott, 1995]. The potential difficulty 
associated with the OGY method is the long transient 
times before the system enters local region where the 
perturbation is effective [Petrovet al., 1994]. In 
conventional control approach, various methods 
including a proportional-Integral (PI) controller 
[Pellegrini and Biardi, 1990] and a modified PI/PID 
controller [Bandyopadhyay, et al., 1997] have been 
employed for controlling and operating the chaotic 
reactors under favourable conditions. However, the 



complex nature of nonlinear dynamical systems 
severely limits the use of conventional linear 
controllers to provide the desired operating 
performance. In advanced control approach; different 
nonlinear model predictive controllers have been 
reported to stabilize the oscillatory and chaotic 
dynamics in continuous chemical and polymerization 
reactors [Lima, 1996, Qammaret al., 1996]. 

Nonlinear control is a class of advanced control in 
which the nonlinear process model serves as the basis 
for controller design. This type of control is expected 
to provide improved control performance since the 
control structure preserves the nonlinearities of the 
real process in the form of its mathematical model. The 
generic model control (GMC) introduced by Lee and 
Sullivan [1988], the globally linearizing control (GLC) 
proposed by Kravaris and Chung [1987] and the 
nonlinear internal model control (NIMC) proposed by 
Henson and Seborg [1991] are the prominent nonlinear 
model based controllers. The GMC allows the 
implementation of the nonlinear process model 
directly into the controller structure. The GLC is a 
model based controller that transforms a nonlinear 
input/output system into a linear input/output system 
through a nonlinear transformation. The NIMC 
proposed by Henson and Seborg is different from the 
GMC of Lee and Sullivan [1988] and the GLC of 
Kravaris and Chung [1987] in that it includes implicit 
integral action in the control structure by using the 
difference between the process output and model 
output as a feedback signal. A nonlinear filter is 
employed in NIMC which provides a tuning 
parameter that can be adjusted for process/model 
mismatch. Since reasonably accurate mathematical 
models have been developed to predict the 
complicated dynamics of nonlinear systems, design of 
nonlinear controllers based on these models is a useful 
control alternative for complex dynamical systems. 
This work is focused towards the model based control 
of nonlinear dynamical systems in the frame work of 
NTMC approach. The NIMC approach based 
controllers have been reported earlier for nonlinear 
systems where the dynamic behaviour is not so 
complex [Kurtz and Henson, 1997, Venkateswarlu and 
Gangiah,1997]. 

The objective of this paper is to derive a controller 
based on NIMC approach for complex dynamical 
systems that incorporates the nonlinear model 
structure and the estimator dynamics in the controller 
formulation. A model based estimator is designed to 
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provide the unmeasured process states that capture 
the fast changing nonlinear dynamics of the process to 
incorporate in the controller. The estimator uses the 
mathematical model of the process in conjunction with 
the known process measurements to estimate the 
states. The design and implementation of the 
proposed estimator based NIMC strategy is studied by 
choosing two typical continuous chemical and 
polymerization reactors that present challenging 
operational and control problems due to their complex 
open-loop dynamics such as input-output 
multiplicities, parametric sensitivity, nonlinear 
oscillations and chaos. Further the comparison of the 
estimator based NIMC strategy is made with 
conventional controllers. 

Control Algorithm 

The general form of a single input single output (SISO) 
system with state space description is: 



x = f(x) + g(x)u 
y = Kx) 



(i) 

(2) 



where x is the vector of states, u is the manipulated 
input, y is the measured output, fix) and g(x) are 
vector functions and h(x) is the scalar function. The 
relative order of the system defined by Eqs. (1 ) and (2) 
is expressed by the following relations: 



L L, h(x) = 0,k<r-l 



L L/~ l h(x) ^ 



(3) 
(4) 



where r represents the relative order of the system and 
Lfh(x) is the Lie derivative of the scalar function h(x) 
with respect to the vector function f(x) with LP/ 
h(x)=h(x). Similarly higher order Lie derivatives as well 
as the Lie derivative of the scalar function h(x) with 
respect to the vector function f(x), and then with 
respect to the vector function g(x) can be defined. The 
relative order defined by Eqs. (3) and (4) represents 
the number of times the output y must be 
differentiated with respect to time so that the input u 
appears explicitly. 

A general form of the control law is written as 
u=p(x) +q(x) v (5) 
where v is the new input. 

The NIMC approach includes an implicit integral 
action by using the difference between the plant and 
model outputs as a feedback signal: 



e = y d -(y-y) (6) 

where y is the process model output. The feedback 
signal is simplified to e=yd with perfect model 
assumption. A nonlinear filter is employed in NIMC 
which provides a tuning parameter that can be 
adjusted for process/model mismatch. The new input v 
for NIMC is defined as 



(7) 



v = -T r L r f l h{x) - t r _ x L r f 2 h(x)-.... 

-T 1 h(x) + r 1 e 
where Ti are controller tuning parameters; r is the 
relative order and e=yd. According to this approach, 
the control law for the system is given by 



v - L f r h(x) 
U = — = p(x) + q(x)v 



L g L/ %(x) 
Application Processes 

Two typical continuous nonlinear dynamical 
processes are chosen for the design and 
implementation of the estimator based NIMC strategy. 

Che mica I re actor 

A non-isothermal, irreversible, first order series 
reaction A— >B— >C in a CSTR with control input and 
load disturbances is described by the following 
dimensionless mass and energy balance equations 
[Kahlertet al.,1981]: 



(8) 



dx 



dt 

dx 



— = 1 - Xy - Da Xy exp(x^ j\ + £ ^x^ ) + dy 



dt 



= + Da Xy exp(xj j\ + £ ^x^ ) 



(9) 



(10) 



Da Sx2exp(fex^ j\ + £ ^3 ) + 



-Xj +B Daxy exp(^ j\ + £ ^3) 



dx^ 
dt 

-DaBa S*2 ex P(^3 A + £ a x 3 ) ~ P^ x ^ ~ x ^ c )(H) 

The variables xi and xi denote the dimensionless 
concentrations of species A and B respectively, and xa 
is the dimensionless temperature. The parameter X3c 
represents the reactor coolant temperature. An 
externally manipulated variable Ut can be defined to 
denote a measure of the deviation in the coolant 
temperature from the reference value X3c. The load 
disturbances in feed compositions are denoted by d\ 
and di, and the load disturbance in the reactor 
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temperature is denoted by da. This reactor system 
exhibit multi-stationary behaviour, oscillations and 
chaos for the parameter values shown in Table 1. 
TABLE 1 CHARACTERIZATION OFSTEADYSTATES 



Set 
No 


Parameter values 




X2s 


X3s 


Stability 


I 


Da=0.06,S=0.0005,eA=0, 
k=l,a=0.426, |3=7.7,B=55.0 


0.0378 


0.9501 


6.0500 


Unstable 


n 


Da=0.26,S=0.5,eA=0,k=l, 
ofO.426, (3=7 7,B=57 .77 


0.0729 


0.1259 


3.8902 


Stable 
limit 
cycle 


ni 


Same as Set II except 
(3=7.9999 


0.0819 


0.1391 


3.7627 


Chaotic 



Homopolymerization reactor 

Continuous homopolymerization reactors represent 
one important class of dynamical systems. The 
homopolymerization reactions are characterized by 
high heat release, high viscosity and poor heat transfer. 
Vinyl acetate homopolymerization has shown 
complex dynamic behaviour for certain ranges of 
operating conditions due to the large heat of 
polymerization, the gel effect and the large activation 
energy of the initiation step. These factors can readily 
produce phenomena such as multiple steady states, 
steady state instability, limit cycles and chaos. 

The kinetic mechanism for free-radical vinyl acetate 
homopolymerization reaction is as follows. 



-> 2R (Radical Initiation) 



R + M- 
P,+M 
P. + P. 



(12) 



-> P i+1 (Chain Propagation) 



->A, 



(Chain Termination) 



The mathematical model for homopolymerization 
reaction of vinyl acetate in tertiary butanol using 
AIBN as initiator is described by the following 
equations [Pinto and Ray, 1995]: 

Monomer Mass Balance 
dv m _ P mf v mf _ (MW) m RJ 



dT dd/pj 



(13) 



+ v p 



dr dT 



Solvent Mass Balance 
dv s PsfV sf 



dr p s 

Initiator Mass Balance 



dT d(l/p s ) 



dr dT 



dc : 



i L = c i f ~ c Mo<ii +k d 0) 

dr 

Energy Balance 

dT = P f (T f -T) _e 

dr " p + pC p (MW) fj 



(15) 



AH 



V mPm k p P 



6 UA(T-T C ) 

pc p V 



(16) 



+ J3u, 



In the above equations, q Q q i is a term that takes into 
account changes in the density of the reactive mixture 
and is denned by 



Mi 



dT 



■1 



(17) 



dr 

and fi = 



dm dm do/ Pf) 

rs s jrr, rm m ,m rpp 



dT 

e ua 



dT 



dT 



pc P V 

The overall density and specific heat of the mixture are 

P = P m v m +P s v s +P p v p 
C - y C + y C + y C 

P J m P m J s P s J P PP 

Under quasi steady state assumption for the free- 
radical species (R and P), the normalized rate of 
monomer reaction is expressed as 



Rate m =k p P + R ckr 



(18) 



where concentration of species P can be written as 

(19) 



The R ckr is the normalized rate of initiation as given 
by 



R ckr ~ 



M 



(20) 



By substituting the value of Rckr from Eq. (20) in Eq. 
(18), we have 
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Rate m =k p P + 



2 f k d c i 



M 

P v 

^ in m 



(21) 



where R m is the rate of monomer consumption 
expressed by 



p Hm m + 

m (MW) m d ' 



(22) 



The termination and propagation constants (k t ,k p ) 
used in the above equations include gel effect with the 



form 

*, = Kg, 
K=Ks P 



(23) 



The gel effect correlations are 

g t = exp(-0.4407x, -6.753* ( 2 -0.3495*, 3 ) 

In Eq. (23), k t and k are the kinetic constants at zero 
polymer concentration, g t and g are the gel effect 
TABLE 2 VA HOMOPOLYMERIZATION SYSTEM PARAMETERS 



f= 0.8 






V = 500 (ml) 




UA = 12(cal/minK) 




AH = 21000(cal/grrol) 




Cpm = 0.470 (cal/gK) 




Cps = 0.716 


(caVg K) 




Pm(T) = 


958.4-1.3276(7- 


273) (g/lt) 


P P (T) = 


1211-0.8496(7- 


273) (g/lt) 


P S {T) = 


74120 ^ 


(g/lt) 


^ (60.21 + 0.1 167) j 




C pp = 0.3453 + 9.55 xl0" 4 (r- 298) (caVgK) 




-6100 




k° p =82.212 xl0 8 e fir (lt/gmol.s) 




-30800 




^ =94.8xl0' 5 e Rr (lVgrrols) 




-2462 




=469.392 xl0 10 e sr (lt/gmoLs) 



correlations that take into account the effects of 
increasing polymer concentration. The other 
parameters necessary for simulation are presented in 
Table 2. 

Design of Estimator Based Nonlinear 
Controller 

The nonlinear internal model controller (NIMC) relies 
on mathematical model of the process incorporating 
process state variable information in the controller 
formulation. Since the state variables desired by the 
controller are not easily available through 
measurement or available with large measurement 
delays, a model based nonlinear estimator is designed 
to provide the unmeasured states needed by the 
controller. 

Controller Des ign 

The design of NIMC for complex dynamic processes 
considered in this work is briefed as follows: 

1) Chemical Reactor 

The fix) and g(x) in Eq. (1) can be written by 
omitting the load disturbances from the 
dimensionless mass and energy balances in Eqs. (9) 
-(11): 



1 - *j - Dax^ exp 



1 + £ A X 3J 



- ^ + Dax^ exp 



- DaSx^ exp 



fct„ 



l + e A x 3 



r -3 ^ 








- BDaaSx^ exp 




1 + £ A x n 






V A 3/ 




V A 3J 



8(x) = 



(24) 



(25) 



The controlled output in Eq. (2) is given by 

y - h(x) - x 3 (26) 

From Eqs. (3) and (4), the relative order of the 
system is evaluated as one. Thus for r =1, the new 
input v in Eq. (7) is obtained as 



V = Tj {X 'l - X 3 ) 



(27) 



On computing the Lie derivatives and substituting 
them along with the new input in Eq. (8) leads to 
the following control law 
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x 3 + BDax x exp 



y l+£ A X 3 j 



BDaaSx 2 exp 



K l + s A x 3j 



/3x 3 



(28) 



where 71 is the NIMC tuning parameter. 
2) Polymerization reactor 



The/fx) and g(x) in Eq. (1) can be defined from the 
mass and energy balances of Eqs. (13) - (16): 



/(*) = 



g(x) 



( 



W Wm 



p 9 9 

m 

P sfsf %Vs 



d(l/p ) 

m m +v p in 



m' m <]j 



P f( T f- T ) 9 AH ,9 UA(T-T) 

v A v p k P - 



pC (MW) m' m p pQ y 
p m p 



dd/p s ) 



p 9 9 s s dT 



if 1 Q- 

9 1 d 1 



P A T f- T ) 9 AH 



pJT -T\ AHv p k P , T . (T-T) 
A / A m m P UA K c' 



pC (MW) V pC 
p m p 



-v p k P 



9 UMT-T) 



pC (MW) m- m p pQ y 
p m p 



(30) 

The controlled output in Eq. (2) is obtained as 
y = h(x) = T (31) 
The conditions in Eqs. (3) and (4) defining the 



(29) 



polymer reactor system is of relative order one. The 
new input v in Eq. (7) is computed as 



v = rAT-T) 



(32) 



Tl (T s -T)- 



p f (T f -T) AHv n 



k p P n 



P UA(T-TJ 



u, = 



Op 



pC p (MW) m V pC p 



(33) 



Estimator Design 

The success of estimator based controller for nonlinear 
dynamical systems relies on the performance of the 
estimator as well as the model of the process that 
supports the estimator. Since mathematical models 
with enough details are used to develop the NIMC 
strategy for complex dynamical systems, it is 
advantageous to use the same models in the estimator 
design in order to obtain fast and accurate estimation 
of unmeasured process states desired by the controller. 
Recently, model based estimators have been 
successfully employed for state estimation in various 
stable systems [Venkateswarlu and Gangiah, 1992, 
Schuler and Schmidt, 1993, Sargantanis and Karim, 



1994, Aguilar-Lopez and Martinez-Guerra, 2005]. In 
this work, a nonlinear model based estimator known 
as extended Kalman filter (EKF) is designed to obtain 
the state estimates from the known temperature 
measurements of the respective processes. The 
estimator design for NIMC strategy of this work is 
similar to that of the recently reported estimator 
involved in the GLC control of a chaotic chemical 
reactor [Karri et al., 2009, Karri, 2011]. The general 
process representation for model based state 
estimation is given in Appendix A. The EKF 
estimation algorithm is given in Appendix B. The 
elements of the state transition and measurement 
matrices involved in the estimator are given in 
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Appendix C . 

Results and Discussion 

Chemical Reactor 

The reactor system described by Eqs. (9)-(ll) exhibits 
multi-stationary behaviour, oscillations and chaos for 
parameter values given as sets I, II and III, respectively, 
in Table. 1. For parameter set I, the system shows 
unstable steady state, where as for parameter sets II 
and III, the system exhibits limit cycle oscillations and 
chaotic behaviour. The temperature data used for state 
estimation are obtained through numerical integration 
of model equations using Gear's method with a 
sampling time of 0.0001 units. A model based state 
estimator, extended Kalman filter is employed to 
estimate the states xi, xi and X3 using the temperature 
measurement of the reactor. The performance of the 
EKF estimator is evaluated by considering different 
cases of process parameter values shown in Table 1. 
The results in Fig. 1 represent the phase plane plots of 
the actual and estimated states of CSTR corresponding 
to the parameter values in set III. 



applied with the objectives of controlling the chemical 
reactor at an unstable steady state in the multiplicity 
region, at the unique unstable steady state responsible 
for limit cycle oscillations and at unique unstable 
steady state responsible for chaotic motion. The 
controller parameter involved in NFMC is tuned as 
n=14.0 based on objective function to minimize the 
integral of absolute error (IAE). The results are also 
compared with those of a modified PID controller 
with feedback mechanism, which has the form given 
by [22], 



du 



dt 



dc c 

L = sJ {t)e + sJ(t) — + sAedt 

dt J 



(34) 



where f(t)=t and f'(t)=l. The tuning parameters of the 
PID controller (a>,&) are evaluated using Ziegler- 
Nichols method. A sampling time of 1 sec is used for 
the implementation of the estimator and the controller. 
For parameter values in set I, the system exhibits 
multiple steady state behaviour. The controller goal is 
to shift the process operating at an arbitrary point (xw 
= 0.04, X2o = 0.9, X3o =5.75) to the set condition in Table 1 
and maintain it at that state. 



0.18 
0.16 
0.14 
0.12 ■ 

0.1 
0.08 
0.06 
0.04 
0.02 




0.18 
0.16 
0.14 
0.12 
0.1 
0.08 
0.06 
0.04 
0.02 





FIG. 1. ACTUAL , X 2 ) AND ESTIMATED 
PHASE PLANE PLOTS FOR SET - m 



X 2 ) 



Fig. 2 shows the actual and estimated concentration 
profiles corresponding to set EI in the presence of 
stochastic noise in temperature obeying Gaussian 
distribution. Similar results are evaluated for set I and 
set II. In all the cases, the estimated states are in close 
resemblance with the actual states. 

The NIMC strategy supported by the state estimator is 








5 6 7 

Time [mm] 




5 6 7 

Time [mm] 



FIG. 2 ACTUAL (jiq ) AND ESTIMATED ( 

CONCENTRATION PROFILES FOR SET-HI IN THE PRESENCE 
OF STOCHASTIC DISTURBANCE 

Fig. 3 shows the process output and controller output 
plots of NIMC and PID controller, respectively for this 
condition. For parameters in set II, the system exhibits 
sustained oscillatory behaviour (stable limit cycle) and 
the controllers are also required to regulate the 
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trajectory at the unique unstable steady state, which is 
the desired condition. Fig. 4 shows the results of the 
NIMC and PID controller to realize the set condition in 
Table 1 from an arbitrarily selected initial state (xio = 
0.08, xio = 0.103, X3o=3.654). For parameters in set IE, the 
system exhibits chaotic behaviour. In this objective, the 
controller has to stabilize the chaotic trajectory exactly 
at the corresponding unique unstable steady state. To 
realize this objective, the controllers are also applied 
using the same initial condition as in set II to satisfy 
the desired condition for set EI given in Table 1. 

The process output and controller output plots of 
NTMC and modified PID controller corresponding to 
this objective are shown in Fig. 5. These plots indicate 
that the controllers fulfil the desired objectives without 
any offset. The nonlinear controller has shown better 
stabilization over PID controller. The performance of 
the controllers is also evaluated by applying them for 
controlling the system at the desired conditions in the 
presence of deterministic and stochastic load 



disturbances. A deterministic disturbance of 1.0 and a 
stochastic load disturbance generated through random 
Gaussian noise of zero mean and a standard deviation 
of 0.25 are considered to represent the ds in 
temperature measurement. The results of the 
controllers for stabilizing the system at unstable steady 
state responsible for chaotic motion with these 
disturbances conditions are shown in Figs. 6 and 7. 

These results indicate the better performance of the 
NIMC strategy in the presence of either type of load 
disturbances. The results of the controllers are also 
evaluated for a series of step changes in the set points 
of the controlled variables. The results in Fig. 8 show 
the process output and controller output plots of 
NIMC and modified PID controller for parameters in 
set IE for a series of ±20% step changes in the set 
conditions in Table 1. These results show the better 
performance of the estimated supported nonlinear 
controller. 




-0.5 -I , , , , , -8 -I , , , , , 

01 234501 23 45 

Time [mill] Time [min] 



FIG. 3 PROCESS AND CONTROLLER OUTPUT PLOTS FORSET I OF (A) NMCAND (B) PID 
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PID 




1.5 
1 

0.5 
- 
-0.5 
-1 

-1 .5 



2 3 
Time [mill] 



- NIMC 



0.5 



-0.5 



-1 .5 




1 2 3 

Time [mill] 



-PID 



2 3 4 

Time [mill] 



FIG. 4 PROCESS AND CONTROLLER OUTPUT PLOTS FOR SET H OF (A) NMC AND (B) PID 



4.5 
4.4 
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Polymerization Reactor 

The mathematical model of homopolymerization 
reactor is solved by numerical integration using Gear's 
method with a sampling time of 0.0001 units. The 
bifurcation diagram of VA homopolymerization 
shown in Fig. 9 drawn for the operating conditions, Tc 
= 318 K, Tf= 315 K, vmf= 0.3 and a; = 0.03203 gmol/lt, 
presents the existence of multiple steady states, limit 
cycle and bifurcation points as well as stable and 
unstable zones with respect to the variation of 
residence time. 

From Fig. 9, it can be observed that the transition 
occurs from stable to unstable zone beginning at 9 = 24 
min. This transition causes oscillatory behaviour in 
response that eventually lead to period doubling 
bifurcation cascade and chaos as shown in Fig. 10(a). 
The magnification results of Fig. 10(a) are further 
explained in Figs. 10(b) and 10(c), which show the 
sequence of periodic windows of periods 2, 4, 8 and so 
on leading to chaotic region. At 9 = 27.37 min, the 
chaotic behaviour with huge temperature oscillations 
is observed as shown in Fig. 11. The time domain 
response is shown in Fig. 11(a) and the corresponding 



phase plane plot is depicted in Fig. 11(b). 
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FIG. 9 PERIOD DOUBLING BIFURCATION DIAGRAM 

A model based estimator, extended Kalman filter (EKF) 
is employed to estimate the states Vm vs a and T using 
the temperature measurement of the reactor. The 
performance of the estimator is evaluated under 
different conditions. Fig. 12 compares the actual and 
estimated results of VA homopolymerization reactor. 
These results confirm the usefulness of the method of 
EKF as a state estimator for polymerization reactor. 
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The NIMC strategy supported by the state estimator is 
applied with the objectives of controlling the 
homopolymerization reactor under different 
conditions and the results are presented in terms of 
three cases. In the first case, the process is initially 
operated under stable operation. Then changes are 
introduced in temperature set point and the controller 
is applied to maintain the desired operation under 
different input disturbance conditions introduced at 
different timings. In the second case, the process is 
operated under oscillatory condition. The controller is 
applied to stabilize the process at the desired 
condition at a time during the open-loop oscillatory 
dynamics. In the third case, the process is operated 
under chaotic behaviour. The controller is applied to 
suppress the chaotic dynamics and to maintain the 
desired operation. The results of the NIMC strategy 
are also compared with those of a conventional PI 
controller. The tuning parameter in NIMC strategy is 
chosen as n =16.5. The PI controller parameters are 
evaluated by using Ziegler and Nichols method and 
further tuned and set as Kc = 1.76 and a =4.0. A 
sampling time of 1 sec is used for implementing the 
estimator and controller. For the first case, the process 
is initially at the stable steady state with T=366.4 K and 
Ti=320 K. At time t = 8.3 min, the temperature set-point 
is changed from 366.4 to 364.87 K and at t = 200 min, Tf 
is changed from 315.00 to 313.00 K while keeping the 
set-point as 364.87 K. The input and output profiles 



corresponding these conditions are shown in Figs. 13(a) 
and 13(b), respectively. In the same case, at t = 400 min, 
the monomer volume fraction, Vmf is changed from 0.3 
to 0.33 while keeping the set-point at 364.87 K and T/at 
313 K. The input and output profiles corresponding 
these conditions are shown in Figs. 13(c) and 13(d), 
respectively. Though the changes are introduced at 
different timings in one time scale, the results are 
represented in separate figures with two time scales 
for the sake of clarity. From these results, it is observed 
that the NIMC tracks the desired condition very 
quickly, whereas the PI controller exhibits oscillatory 
behaviour for a longer duration. The PI controller 
could even fail for change in Vmf as shown in Fig. 13(c). 
The NIMC effectively rejects the disturbance 
conditions and leads the process back to its desired 
operation in a very conservative fashion. 

In the second case, the process is disturbed from its 
initial steady state by changing the Tc from 320 K to 
318 K at t=8.3 min. This condition makes the process to 
exhibit oscillatory behaviour. This open-loop 
behaviour of the system is shown in Fig. 14(a). LXiring 
the oscillatory operation at f = 217.4 min, the 
controllers are applied to stabilize the operation and to 
maintain the desired response condition, which is set 
as 364.87K for this case. The results of process output 
and controller output for both the controllers are 
shown in Figs. 14(b) and 14(c), respectively. From 
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these results, it is observed that both the controllers 
are able to track the process to the desired condition, 
but the NIMC has shown better performance by 
quickly suppressing the oscillations. In the third case, 
the process is exhibiting chaotic behaviour for T=347.6 
KandT r = 318K. 
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FIG. 14 OPEN LOOP AND CLOSED LOOP RESULTS: (A) OPEN- 
LOOP BEHAVDUR (B) CONTROLLEDOUTPUT(C) CONTROL 
INPUT 

During the chaotic operation at t = 240 min, the 
controllers are applied to stabilize the operation and to 
maintain the desired response which is set as 364.87K 
for this case. The results of process output and 
controller output for both the controllers are shown in 
Figs. 15(a) and 15(b), respectively. Again, the NIMC 
has performed better by quickly suppressing the 
chaotic dynamics with less stringent control actions. 

The results evaluated for both the chemical and 
polymerization reactors show the superior 
performance of the estimator based NIMC strategy 
over the conventional controllers. The performances of 
the controllers for different conditions are also 
expressed in terms of IAE values. The IAE values 



shown in Table 3 correspond to the conditions for 
which the results are depicted in figures referred in 
Table. 3. These quantification results further confirm 
that the NTMC strategy is better suited for the control 
of nonlinear dynamical processes. 
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FIG. 15 OPEN LOOP AND CLOSED LOOP RESULTS : (A) 
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TABLE 3 CONTROLLER'S PERFORMANCE 



Chemical Reactor 


IAE 


REFERENCE 


NIMC 


PID 


0.938 


365.163 


Fig. 3 


0.748 


72.763 


Fig. 4 


0.347 


93.524 


Fig. 5 


7.121 


83.042 


Fig. 6 


0.924 


93.709 


Fig. 7 


6.852 


118.350 


Fig. 8 


Polymer Reactor 


0.000985 


174.193 


Fig. 13 


40.605 


110.697 


Fig. 14 


138.626 


339.425 


Fig. 15 



Conclusions 

Complex dynamical systems can present the 
phenomena such as multiple steady state behaviour, 
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nonlinear oscillations and chaos. This phenomenon is 
usually analysed by means of dynamic mathematical 
models. The availability of such models can form a 
basis for the development of advanced nonlinear 
controllers for such systems. The nonlinear model 
based controllers require the process state information 
in the controller formulation. This necessitates the 
need for an estimator to provide the unmeasured 
process states desired by the controller. Thus, this 
work is focused towards the development of a 
nonlinear model based controller for complex 
dynamic systems, which is supported by a nonlinear 
model based estimator. A nonlinear internal model 
control (NIMC) strategy that incorporates the 
nonlinear model structure and the estimator dynamics 
in the control law is presented for the control of 
complex dynamical systems. The estimator is designed 
to provide the unmeasured process states that capture 
the fast changing nonlinear dynamics of the process to 
incorporate in the controller. The performance of the 
estimator supported NIMC strategy is evaluated by 
applying it for the control of a non-isothermal 
nonlinear chemical reactor and a homopolymerization 
reactor, which exhibit rich dynamical behaviour 
ranging from stable situations to chaos. The results 
evaluated under different conditions show the 
superior performance of the estimator based NIMC 
strategy over the conventional controllers. 

Appendix A: Process Representation 

The time varying model of the nonlinear process is 
represented by 



x(t) = f(x(t),t) + w(t) , x(0) = x 



(A.1) 



where x(t) is n dimensional state vector, / is a 
nonlinear function of state x(t) and w(t) is an additive 
Gaussian noise with zero mean. The nonlinear 
measurement model with observation noise can be 
expressed as 



y(t k ) = h(x(t k )) + v(t k ) 



(A.2) 



where h is a nonlinear function of state x(tk). The 
expected values of noise covariance matrices for the 
initial state x(0), process noise w(t) and observation 
noise v(tk) are given by the following relations, 



P =E[(x -x(0))(x -x(0)) T ] 
£(0 = £[w(t)w T (t)] 
i?(f,) = £[v(t k )v T (t k )] 



(A.3) 



where Pois initial statecovariance matrix, Q (t) is 
process noise covariance matrix and R(tk) is 
observation noise covariance matrix. The EKF 
estimation algorithm is given in Appendix B. 

Appendix B: State Estimation Algorithm 

The extended Kalman filter is computed in two steps. 
The first is a prediction step, which is used to 
extrapolate the previous best estimates, and the 
second is a correction step by which the updated 
estimates are formed. Since prediction is based on 
process model, continuous prediction and discrete 
correction are employed in the estimation scheme. 

Prediction equations: By starting with an initial estimate 
Xo and its covariance Po at time zero and no 
measurements are taken between tk-i and fe , the 
propagating expression for the state estimate and its 
covariance from fjt-i to fe are, 



x(t/t k _ l ) = f(x(t/t k _ l ),t) 

P(t/t k _ 1 ) = F(x(t/t kA ),t)P(t/t k _ l ) 

+ P(t/t kA )F T (x(t/t kA ),t) + Q(t) 

where F(x(t / t k _ l ),t) is the state transition matrix 
whose ij th element is given by 

df,(x(t),t) 



(B.l) 
(B.2) 



F{x{t/t kA ),t) = 



d Xj (t) 



x(t)=x(t/t lA ) 



(B.3) 



The solution of the propagated estimate x{t/t kA ) and 
its covariance P(t/tk-i) at time tk are denoted by 
x(t k /t kA ) and P(tk/tk-i). By using measurements at time 

tk, the update estimate x(t k /t k ) and its covariance 
P(tk/tk) are computed. 

Correction equations: The equations to obtain corrected 
estimates are: 

x(t k /t k ) = x(t k /t k _ x ) + K(t k )[y(t k ) - h (x(t k I t k _ x ))] 

(B.4) 

P(t k /t k ) = (/ - K(t k ) H(x(t k )))P(t k I t k _ x ) 



P{t k lt k _,)H T (x{t k )) 



(H(x(t k ))P(t k /t k _ 1 )H T (x(t k )) + R) 
dh t (x(t k y) 



(B.5) 
(B.6) 



where, H(x(t k )) = 



dx(t k ) 



The recursive initial conditions for state and 
covariance are: 



15 



www.se ipub .org/acer 



Advanced Chemical Engineering Research (ACER) Volume 2 Issue 1, March 2013 



P(t k /t k _ l ) = p(t k /t k ) 



(B.7) 



Appendix C: Elements of F and H Matrices Chemical 
reactor 

The nonlinear dynamic model of the reactor system in 
its dimensionless form is used in conjunction with the 
temperature measurements to estimate the reactor 
species concentrations. The design of soft sensor 
involves the following components: 

State vector 

The chemical species concentrations and the 
temperature define the state vector as 



X — L-^'j "^"2 "^"3 -J 

State transition matrix 



(C.l) 



fn 


fn 


fn 


fn 


fl2 


fn 


fn 


fyi 


fl3 



The elements fij of the state transition matrix, Fare 
computed by taking the partial derivatives of fix) 
defined by (24) - (26) with respect to the state vector: 



F = 



(C.2) 

Measurement matrix 

The measurement relation for temperature is 

H=[0 1] (C3) 

The temperature state equation, (26) in its discrete 
form is used as the nonlinear measurement equation, 
h(x). The elements of the measurement transition 
matrix, Hx are computed by taking the partial 
derivatives of h(x) with respect to the state vector 



H x =[h n h n h l3 ] 



(C.4) 



All these components are evaluated for chaotic reactor 
and used with the EKF estimator to obtain measured 
and unmeasured states of the reactor. The soft sensor 
uses the temperature data of every sampling time as 
measurements and provides the estimates of 
temperature as well as reactor species concentrations. 

Polymer reactor 

The nonlinear dynamic model of the VA 
homopolymerization system is used in conjunction 
with the temperature measurements to estimate the 
volume fraction of monomer, solvent, concentration of 



initiator and reactor temperature. The design of soft 
sensor involves the following components: 

State vector 

The volume fraction of monomer, solvent, 
concentration of initiator and reactor temperature 
define the state vector as 



x = lv m v s c t T] 
State transition matrix 



(C.5) 



The elements fj of the state transition matrix, Fare 
computed by taking the partial derivatives of fix) 
defined by (28) - (30) with respect to the state vector: 



fn 


fn 


fn 


fn 


fn 


fn 


/* 


/* 


fn 


fx 


fa 


/ 3 4 


/« 


f\2 


U 


fw 



F = 



Measurement matrix 

The measurement relation for temperature is 
H =[0 1] 



(C.6) 



(C.7) 



The temperature state equation, (30) in its discrete 
form is used as the nonlinear measurement equation, 
h(x). The elements of the measurement transition 
matrix, Hx are computed by taking the partial 
derivatives of h(x) with respect to the state vector 



\h n h l2 h n 



(C.8) 



All these components are evaluated for 
homopolymerization reactor and used with the EKF 
estimator to obtain measured and unmeasured states 
of the reactor. The state estimator uses the temperature 
data of every sampling time as measurements and 
provides the estimates of temperature as well as 
reactor species concentrations. 

Notation 

B heat of reaction parameter 

a initiator concentration in reactor, gmolll 

af initiator feed concentration, gmolll 

Cp heat capacity of reaction mixture, cal/g °C 

c P f heat capacity of feed stream, cal/g °C 

di, di, di load disturbance in xi, X2, X3 

Da Damkohler number 

e set point error 
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/ initiator dissociation efficiency 

-AHr heat of polymerization, cal/gmol 

k ratio of the activation energies for the series 
reaction 

kd rate constant for initiator decomposition, lis 

kp rate constant for propagation, llgmol s 

h rate constant for termination, llgmol s 

Kc controller gain 

(MW)m molecular weight of monomer, glgmol 

P live radical concentration, gmolll 

qoCji ratio of outlet to inlet volumetric flow rate 

Km rate of consumption of monomer, gmolll min 

S ratio of the rate constants for the series 

reaction 

t time, min 

T reactor temperature, K 

Tc coolant temperature, K 

Tf feed temperature, K 

UA overall heat transfer coefficient, callmin K 

Ut manipulated variable 

V volume of the reactor, I 

xt mole fraction of polymer 

xi dimensionless concentration of species A 

X2 dimensionless concentration of species B 

X3 dimensionless temperature 

X3c dimensionless coolant temperature 

xi xi , X2 xt , xs^'set point for variable xi, X2 and X3 

xis , X2s , X3ssteady state values of xi, X2 and X3 

y; weight fraction of respective components 

(i=m,s,p) 

Greek symbols 

a ratio of heat effects for the series reaction 

P heat transfer coefficient 

eo initial controller gain 

ei initial controller gain 

a activation number 

residence time, min 



t dimensionless time (t/0) 

id time constant for derivative action 

v time constant for integral action 

Vm volume fraction of monomer in reactor 

Vmf volume fraction of monomer in feed 

p density of reaction mixture, gll 

PS density of feed, gll 

PUT) density of feed adjusted to reactor 
temperature, gll 

p m density of monomer at reactor conditions, gll 

pmf density of monomer at feed conditions, gll 

Subscripts 

m monomer 

s solvent 

p polymer 

/ feed 
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